home *** CD-ROM | disk | FTP | other *** search
/ CD ROM Paradise Collection 4 / CD ROM Paradise Collection 4 1995 Nov.iso / program / nrpas13.zip / QROMB.PAS < prev    next >
Pascal/Delphi Source File  |  1991-04-29  |  885b  |  35 lines
  1. PROCEDURE qromb(a,b: real; VAR ss: real);
  2. (* Programs using routine QROMB must define type
  3. TYPE
  4.    glnarray = ARRAY [1..n] OF real;
  5. just as for routine POLINT which it calls. In this case
  6. n should have a value no smaller than the constant k below. *)
  7. LABEL 99;
  8. CONST
  9.    eps=1.0e-6;
  10.    jmax=20;
  11.    jmaxp=21;   (* jmax+1 *)
  12.    k=5;
  13. VAR
  14.    i,j: integer;
  15.    dss: real;
  16.    h,s: ARRAY[1..jmaxp] OF real;
  17.    c,d: glnarray;
  18. BEGIN
  19.    h[1] := 1.0;
  20.    FOR j := 1 TO jmax DO BEGIN
  21.       trapzd(a,b,s[j],j);
  22.       IF  (j >= k)  THEN BEGIN
  23.          FOR i := 1 TO k DO BEGIN
  24.             c[i] := h[j-k+i];
  25.             d[i] := s[j-k+i]
  26.          END;
  27.          polint(c,d,k,0.0,ss,dss);
  28.          IF (abs(dss) < eps*abs(ss)) THEN GOTO 99
  29.       END;
  30.       s[j+1] := s[j];
  31.       h[j+1] := 0.25*h[j]
  32.    END;
  33.    writeln('pause in QROMB - too many steps'); readln;
  34. 99:   END;
  35.